home *** CD-ROM | disk | FTP | other *** search
/ Whiteline: delta / whiteline CD Series - delta.iso / tools / anwender / astrolog / sources / charts3.c < prev    next >
C/C++ Source or Header  |  1995-11-25  |  26KB  |  705 lines
  1. /*
  2. ** Astrolog (Version 4.40) File: charts3.c
  3. **
  4. ** IMPORTANT NOTICE: The graphics database and chart display routines
  5. ** used in this program are Copyright (C) 1991-1995 by Walter D. Pullen
  6. ** (astara@u.washington.edu). Permission is granted to freely use and
  7. ** distribute these routines provided one doesn't sell, restrict, or
  8. ** profit from them in any way. Modification is allowed provided these
  9. ** notices remain with any altered or edited versions of the program.
  10. **
  11. ** The main planetary calculation routines used in this program have
  12. ** been Copyrighted and the core of this program is basically a
  13. ** conversion to C of the routines created by James Neely as listed in
  14. ** Michael Erlewine's 'Manual of Computer Programming for Astrologers',
  15. ** available from Matrix Software. The copyright gives us permission to
  16. ** use the routines for personal use but not to sell them or profit from
  17. ** them in any way.
  18. **
  19. ** The PostScript code within the core graphics routines are programmed
  20. ** and Copyright (C) 1992-1993 by Brian D. Willoughby
  21. ** (brianw@sounds.wa.com). Conditions are identical to those above.
  22. **
  23. ** The extended accurate ephemeris databases and formulas are from the
  24. ** calculation routines in the program "Placalc" and are programmed and
  25. ** Copyright (C) 1989,1991,1993 by Astrodienst AG and Alois Treindl
  26. ** (alois@azur.ch). The use of that source code is subject to
  27. ** regulations made by Astrodienst Zurich, and the code is not in the
  28. ** public domain. This copyright notice must not be changed or removed
  29. ** by any user of this program.
  30. **
  31. ** Initial programming 8/28,30, 9/10,13,16,20,23, 10/3,6,7, 11/7,10,21/1991.
  32. ** X Window graphics initially programmed 10/23-29/1991.
  33. ** PostScript graphics initially programmed 11/29-30/1992.
  34. ** Last code change made 1/29/1995.
  35. */
  36.  
  37. #include "astrolog.h"
  38.  
  39.  
  40. /*
  41. ******************************************************************************
  42. ** Multiple Chart Scanning Routines.
  43. ******************************************************************************
  44. */
  45.  
  46. /* Search through a day, and print out the times of exact aspects among the  */
  47. /* planets during that day, as specified with the -d switch, as well as the  */
  48. /* times when a planet changes sign or direction. To do this, we cast charts */
  49. /* for the beginning and end of the day, or a part of a day, and do a linear */
  50. /* equation check to see if anything exciting happens during the interval.   */
  51. /* (This is probably the single most complicated procedure in the program.)  */
  52.  
  53. void ChartInDaySearch(fProg)
  54. bool fProg;
  55. {
  56.   char sz[cchSzDef];
  57.   int source[MAXINDAY], aspect[MAXINDAY], dest[MAXINDAY],
  58.     sign1[MAXINDAY], sign2[MAXINDAY], D1, D2, occurcount, division, div,
  59.     fYear, yea0, yea1, yea2, i, j, k, l, s1, s2;
  60.   real time[MAXINDAY], divsiz, d1, d2, e1, e2, f1, f2, g;
  61.  
  62.   /* If parameter 'fProg' is set, look for changes in a progressed chart. */
  63.  
  64.   fYear = us.fInDayMonth && (Mon2 == 0);
  65.   division = (fYear || fProg) ? 1 : us.nDivision;
  66.   divsiz = 24.0 / (real)division*60.0;
  67.  
  68.   /* If -dY in effect, then search through a range of years. */
  69.  
  70.   yea1 = fProg ? Yea2 : Yea;
  71.   yea2 = fYear ? (yea1 + us.nEphemYears - 1) : yea1;
  72.   for (yea0 = yea1; yea0 <= yea2; yea0++) {
  73.  
  74.   /* If -dm in effect, then search through the whole month, day by day. */
  75.  
  76.   if (us.fInDayMonth) {
  77.     D1 = 1;
  78.     if (fYear) {
  79.       Mon2 = 1; D2 = DayInYearHi(yea0);
  80.     } else
  81.       D2 = DayInMonth(fProg ? Mon2 : Mon, yea0);
  82.   } else
  83.     D1 = D2 = Day;
  84.  
  85.   /* Start searching the day or days in question for exciting stuff. */
  86.  
  87.   for (Day2 = D1; Day2 <= D2; Day2 = AddDay(Mon, Day2, yea0, 1)) {
  88.     occurcount = 0;
  89.  
  90.     /* Cast chart for beginning of day and store it for future use. */
  91.  
  92.     SetCI(ciCore, fYear ? Mon2 : Mon, Day2, yea0, 0.0, Dst, Zon, Lon, Lat);
  93.     if (us.fProgress = fProg) {
  94.       is.JDp = MdytszToJulian(Mon2, DD, yea0, 0.0, Dst, Zon);
  95.       ciCore = ciMain;
  96.     }
  97.     CastChart(fTrue);
  98.     for (i = 1; i <= cSign; i++) {
  99.       cp2.cusp[i] = house[i];
  100.       cp2.house[i] = inhouse[i];
  101.     }
  102.     for (i = 1; i <= cObj; i++) {
  103.       cp2.obj[i] = planet[i];
  104.       cp2.dir[i] = ret[i];
  105.     }
  106.  
  107.     /* Now divide the day into segments and search each segment in turn. */
  108.     /* More segments is slower, but has slightly better time accuracy.   */
  109.  
  110.     for (div = 1; div <= division; div++) {
  111.  
  112.       /* Cast the chart for the ending time of the present segment. The   */
  113.       /* beginning time chart is copied from the previous end time chart. */
  114.  
  115.       SetCI(ciCore, fYear ? Mon2 : Mon, Day2, yea0,
  116.         DegToDec(24.0*(real)div/(real)division), Dst, Zon, Lon, Lat);
  117.       if (fProg) {
  118.         is.JDp = MdytszToJulian(Mon2, DD+1, yea0, 0.0, Dst, Zon);
  119.         ciCore = ciMain;
  120.       }
  121.       CastChart(fTrue);
  122.       for (i = 1; i <= cSign; i++) {
  123.         cp1.cusp[i] = cp2.cusp[i]; cp1.house[i] = cp2.house[i];
  124.         cp2.cusp[i] = house[i];  cp2.house[i] = inhouse[i];
  125.       }
  126.       for (i = 1; i <= cObj; i++) {
  127.         cp1.obj[i] = cp2.obj[i]; cp1.dir[i] = cp2.dir[i];
  128.         cp2.obj[i] = planet[i];  cp2.dir[i] = ret[i];
  129.       }
  130.  
  131.       /* Now search through the present segment for anything exciting. */
  132.  
  133.       for (i = 1; i <= cObj; i++) if (!ignore[i] && (fProg || FThing(i))) {
  134.         s1 = SFromZ(cp1.obj[i])-1;
  135.         s2 = SFromZ(cp2.obj[i])-1;
  136.  
  137.         /* Does the current planet change into the next or previous sign? */
  138.  
  139.         if (!ignore[i] && s1 != s2 && !ignore[0]) {
  140.           source[occurcount] = i;
  141.           aspect[occurcount] = aSig;
  142.           dest[occurcount] = s2+1;
  143.           time[occurcount] = MinDistance(cp1.obj[i],
  144.             (real)(cp1.dir[i] >= 0.0 ? s2 : s1) * 30.0) /
  145.             MinDistance(cp1.obj[i], cp2.obj[i])*divsiz + (real)(div-1)*divsiz;
  146.           sign1[occurcount] = sign2[occurcount] = s1+1;
  147.           occurcount++;
  148.         }
  149.  
  150.         /* Does the current planet go retrograde or direct? */
  151.  
  152.         if (!ignore[i] && (cp1.dir[i] < 0.0) != (cp2.dir[i] < 0.0) &&
  153.           !ignore2[0]) {
  154.           source[occurcount] = i;
  155.           aspect[occurcount] = aDir;
  156.           dest[occurcount] = cp2.dir[i] < 0.0;
  157.           time[occurcount] = RAbs(cp1.dir[i])/(RAbs(cp1.dir[i])+
  158.             RAbs(cp2.dir[i]))*divsiz + (real)(div-1)*divsiz;
  159.           sign1[occurcount] = sign2[occurcount] = s1+1;
  160.           occurcount++;
  161.         }
  162.  
  163.         /* Now search for anything making an aspect to the current planet. */
  164.  
  165.         for (j = i+1; j <= cObj; j++) if (!ignore[j] && (fProg || FThing(j)))
  166.           for (k = 1; k <= us.nAsp; k++) if (FAcceptAspect(i, k, j)) {
  167.             d1 = cp1.obj[i]; d2 = cp2.obj[i];
  168.             e1 = cp1.obj[j]; e2 = cp2.obj[j];
  169.             if (MinDistance(d1, d2) < MinDistance(e1, e2)) {
  170.               SwapR(&d1, &e1);
  171.               SwapR(&d2, &e2);
  172.             }
  173.  
  174.             /* We are searching each aspect in turn. Let's subtract the  */
  175.             /* size of the aspect from the angular difference, so we can */
  176.             /* then treat it like a conjunction.                         */
  177.  
  178.             if (MinDistance(e1, Mod(d1-aspectangle[k])) <
  179.                 MinDistance(e2, Mod(d2+aspectangle[k]))) {
  180.               e1 = Mod(e1+aspectangle[k]);
  181.               e2 = Mod(e2+aspectangle[k]);
  182.             } else {
  183.               e1 = Mod(e1-aspectangle[k]);
  184.               e2 = Mod(e2-aspectangle[k]);
  185.             }
  186.  
  187.             /* Check to see if the aspect actually occurs during our    */
  188.             /* segment, making sure we take into account if one or both */
  189.             /* planets are retrograde or if they cross the Aries point. */
  190.  
  191.             f1 = e1-d1;
  192.             if (RAbs(f1) > rDegHalf)
  193.               f1 -= RSgn(f1)*rDegMax;
  194.             f2 = e2-d2;
  195.             if (RAbs(f2) > rDegHalf)
  196.               f2 -= RSgn(f2)*rDegMax;
  197.             if (MinDistance(Midpoint(d1, d2), Midpoint(e1, e2)) < rDegQuad &&
  198.               RSgn(f1) != RSgn(f2)) {
  199.               source[occurcount] = i;
  200.               aspect[occurcount] = k;
  201.               dest[occurcount] = j;
  202.  
  203.               /* Horray! The aspect occurs sometime during the interval.   */
  204.               /* Now we just have to solve an equation in two variables to */
  205.               /* find out where the "lines" cross, i.e. the aspect's time. */
  206.  
  207.               f1 = d2-d1;
  208.               if (RAbs(f1) > rDegHalf)
  209.                 f1 -= RSgn(f1)*rDegMax;
  210.               f2 = e2-e1;
  211.               if (RAbs(f2) > rDegHalf)
  212.                 f2 -= RSgn(f2)*rDegMax;
  213.               g = (RAbs(d1-e1) > rDegHalf ?
  214.                 (d1-e1)-RSgn(d1-e1)*rDegMax : d1-e1)/(f2-f1);
  215.               time[occurcount] = g*divsiz + (real)(div-1)*divsiz;
  216.               sign1[occurcount] = (int)(Mod(cp1.obj[i]+
  217.                 RSgn(cp2.obj[i]-cp1.obj[i])*
  218.                 (RAbs(cp2.obj[i]-cp1.obj[i]) > rDegHalf ? -1 : 1)*
  219.                 RAbs(g)*MinDistance(cp1.obj[i], cp2.obj[i]))/30.0)+1;
  220.               sign2[occurcount] = (int)(Mod(cp1.obj[j]+
  221.                 RSgn(cp2.obj[j]-cp1.obj[j])*
  222.                 (RAbs(cp2.obj[j]-cp1.obj[j]) > rDegHalf ? -1 : 1)*
  223.                 RAbs(g)*MinDistance(cp1.obj[j], cp2.obj[j]))/30.0)+1;
  224.               occurcount++;
  225.             }
  226.           }
  227.       }
  228.     }
  229.  
  230.     /* After all the aspects, etc, in the day have been located, sort   */
  231.     /* them by time at which they occur, so we can print them in order. */
  232.  
  233.     for (i = 1; i < occurcount; i++) {
  234.       j = i-1;
  235.       while (j >= 0 && time[j] > time[j+1]) {
  236.         SwapN(source[j], source[j+1]);
  237.         SwapN(aspect[j], aspect[j+1]);
  238.         SwapN(dest[j], dest[j+1]);
  239.         SwapR(&time[j], &time[j+1]);
  240.         SwapN(sign1[j], sign1[j+1]); SwapN(sign2[j], sign2[j+1]);
  241.         j--;
  242.       }
  243.     }
  244.  
  245.     /* Finally, loop through and display each aspect and when it occurs. */
  246.  
  247.     for (i = 0; i < occurcount; i++) {
  248.       s1 = (int)time[i]/60;
  249.       s2 = (int)time[i]-s1*60;
  250.       j = Day2;
  251.       if (fYear || fProg) {
  252.         l = Mon2;
  253.         while (j > (k = DayInMonth(l, yea0))) {
  254.           j -= k;
  255.           l++;
  256.         }
  257.       }
  258.       SetCI(ciSave, fYear || fProg ? l : Mon, j, yea0,
  259.         DegToDec(time[i] / 60.0), Dst, Zon, Lon, Lat);
  260.       k = DayOfWeek(fYear || fProg ? l : Mon, j, yea0);
  261.       AnsiColor(kRainbowA[k + 1]);
  262.       sprintf(sz, "(%c%c%c) ", chDay3(k)); PrintSz(sz);
  263.       AnsiColor(kDefault);
  264.       sprintf(sz, "%s %s ",
  265.         SzDate(fYear || fProg ? l : Mon, j, yea0, fFalse),
  266.         SzTime(s1, s2)); PrintSz(sz);
  267.       PrintAspect(source[i], sign1[i],
  268.         (int)RSgn(cp1.dir[source[i]])+(int)RSgn(cp2.dir[source[i]]),
  269.         aspect[i], dest[i], sign2[i],
  270.         (int)RSgn(cp1.dir[dest[i]])+(int)RSgn(cp2.dir[dest[i]]),
  271.         (char)(fProg ? 'e' : 'd'));
  272.       PrintInDay(source[i], aspect[i], dest[i]);
  273.     }
  274.   }
  275.   }
  276.  
  277.   /* Recompute original chart placements as we've overwritten them. */
  278.  
  279.   ciCore = ciMain;
  280.   CastChart(fTrue);
  281. }
  282.  
  283.  
  284. /* Search through a month, year, or years, and print out the times of exact */
  285. /* transits where planets in the time frame make aspect to the planets in   */
  286. /* some other chart, as specified with the -t switch. To do this, we cast   */
  287. /* charts for the start and end of each month, or within a month, and do an */
  288. /* equation check for aspects to the other base chart during the interval.  */
  289.  
  290. void ChartTransitSearch(fProg)
  291. bool fProg;
  292. {
  293.   real planet3[objMax], house3[cSign+1], ret3[objMax], time[MAXINDAY];
  294.   char sz[cchSzDef];
  295.   int source[MAXINDAY], aspect[MAXINDAY], dest[MAXINDAY], sign[MAXINDAY],
  296.     isret[MAXINDAY], M1, M2, Y1, Y2, occurcount, div, i, j, k, s1, s2, s3;
  297.   real divsiz, daysiz, d, e1, e2, f1, f2;
  298.  
  299.   for (i = 1; i <= cSign; i++)
  300.     house3[i] = house[i];
  301.   for (i = 1; i <= cObj; i++) {
  302.     planet3[i] = planet[i];
  303.     ret3[i] = ret[i];
  304.   }
  305.  
  306.   /* Hacks: Searching month number zero means to search the whole year     */
  307.   /* instead, month by month. Searching a negative month means to search   */
  308.   /* multiple years, with the span off the year stored in the "day" field. */
  309.  
  310.   Y1 = Y2 = Yea2;
  311.   M1 = M2 = Mon2;
  312.   if (Mon2 < 1) {
  313.     M1 = 1; M2 = 12;
  314.     if (Mon2 < 0) {
  315.       if (Day2 < 1) {
  316.         Y1 = Yea2 + Day2 + 1; Y2 = Yea2;
  317.       } else {
  318.         Y1 = Yea2; Y2 = Yea2 + Day2 - 1;
  319.       }
  320.     }
  321.   }
  322.  
  323.   /* Start searching the year or years in question for any transits. */
  324.  
  325.   for (Yea2 = Y1; Yea2 <= Y2; Yea2++)
  326.  
  327.   /* Start searching the month or months in question for any transits. */
  328.  
  329.   for (Mon2 = M1; Mon2 <= M2; Mon2++) {
  330.     daysiz = (real)DayInMonth(Mon2, Yea2)*24.0*60.0;
  331.     divsiz = daysiz / (real)us.nDivision;
  332.  
  333.     /* Cast chart for beginning of month and store it for future use. */
  334.  
  335.     SetCI(ciCore, Mon2, 1, Yea2, 0.0, Dst2, Zon2, Lon2, Lat2);
  336.     if (us.fProgress = fProg) {
  337.       is.JDp = MdytszToJulian(MM, DD, YY, 0.0, Dst2, Zon2);
  338.       ciCore = ciMain;
  339.     }
  340.     for (i = 1; i <= oNorm; i++)
  341.       SwapN(ignore[i], ignore2[i]);
  342.     CastChart(fTrue);
  343.     for (i = 1; i <= oNorm; i++)
  344.       SwapN(ignore[i], ignore2[i]);
  345.     for (i = 1; i <= cSign; i++)
  346.       cp2.cusp[i] = house[i];
  347.     for (i = 1; i <= oCore; i++) {
  348.       cp2.obj[i] = planet[i];
  349.       cp2.dir[i] = ret[i];
  350.     }
  351.  
  352.     /* Divide our month into segments and then search each segment in turn. */
  353.  
  354.     for (div = 1; div <= us.nDivision; div++) {
  355.       occurcount = 0;
  356.  
  357.       /* Cast the chart for the ending time of the present segment, and */
  358.       /* copy the start time chart from the previous end time chart.    */
  359.  
  360.       d = 1.0 + (daysiz/24.0/60.0)*(real)div/(real)us.nDivision;
  361.       SetCI(ciCore, Mon2, (int)d, Yea2,
  362.         DegToDec(RFract(d)*24.0), Dst2, Zon2, Lon2, Lat2);
  363.       if (fProg) {
  364.         is.JDp = MdytszToJulian(MM, DD, YY, 0.0, Dst2, Zon2);
  365.         ciCore = ciMain;
  366.       }
  367.       for (i = 1; i <= oNorm; i++)
  368.         SwapN(ignore[i], ignore2[i]);
  369.       CastChart(fTrue);
  370.       for (i = 1; i <= oNorm; i++)
  371.         SwapN(ignore[i], ignore2[i]);
  372.       for (i = 1; i <= cSign; i++) {
  373.         cp1.cusp[i] = cp2.cusp[i]; cp2.cusp[i] = house[i];
  374.       }
  375.       for (i = 1; i <= oCore; i++) {
  376.         cp1.obj[i] = cp2.obj[i]; cp1.dir[i] = cp2.dir[i];
  377.         cp2.obj[i] = planet[i];  cp2.dir[i] = ret[i];
  378.       }
  379.  
  380.       /* Now search through the present segment for any transits. Note that */
  381.       /* stars can be transited, but they can't make transits themselves.   */
  382.  
  383.       for (i = 1; i <= cObj; i++) if (!ignore[i])
  384.         for (j = 1; j <= oNorm; j++) if (!ignore2[j] && (fProg || FThing(j)))
  385.  
  386.           /* Between each pair of planets, check if they make any aspects. */
  387.  
  388.           for (k = 1; k <= us.nAsp; k++) if (FAcceptAspect(i, k, j)) {
  389.             d = planet3[i]; e1 = cp1.obj[j]; e2 = cp2.obj[j];
  390.             if (MinDistance(e1, Mod(d-aspectangle[k])) <
  391.                 MinDistance(e2, Mod(d+aspectangle[k]))) {
  392.               e1 = Mod(e1+aspectangle[k]);
  393.               e2 = Mod(e2+aspectangle[k]);
  394.             } else {
  395.               e1 = Mod(e1-aspectangle[k]);
  396.               e2 = Mod(e2-aspectangle[k]);
  397.             }
  398.  
  399.             /* Check to see if the present aspect actually occurs during the */
  400.             /* segment, making sure we check any Aries point crossings.      */
  401.  
  402.             f1 = e1-d;
  403.             if (RAbs(f1) > rDegHalf)
  404.               f1 -= RSgn(f1)*rDegMax;
  405.             f2 = e2-d;
  406.             if (RAbs(f2) > rDegHalf)
  407.               f2 -= RSgn(f2)*rDegMax;
  408.             if (MinDistance(d, Midpoint(e1, e2)) < rDegQuad &&
  409.               RSgn(f1) != RSgn(f2) && occurcount < MAXINDAY) {
  410.  
  411.               /* Ok, we have found a transit. Now determine the time */
  412.               /* and save this transit in our list to be printed.    */
  413.  
  414.               source[occurcount] = j;
  415.               aspect[occurcount] = k;
  416.               dest[occurcount] = i;
  417.               time[occurcount] = RAbs(f1)/(RAbs(f1)+RAbs(f2))*divsiz +
  418.                 (real)(div-1)*divsiz;
  419.               sign[occurcount] = (int)(Mod(
  420.                 MinDistance(cp1.obj[j], Mod(d-aspectangle[k])) <
  421.                 MinDistance(cp2.obj[j], Mod(d+aspectangle[k])) ?
  422.                 d-aspectangle[k] : d+aspectangle[k])/30.0)+1;
  423.               isret[occurcount] = (int)RSgn(cp1.dir[j]) +
  424.                 (int)RSgn(cp2.dir[j]);
  425.               occurcount++;
  426.             }
  427.           }
  428.  
  429.       /* After all transits located, sort them by time at which they occur. */
  430.  
  431.       for (i = 1; i < occurcount; i++) {
  432.         j = i-1;
  433.         while (j >= 0 && time[j] > time[j+1]) {
  434.           SwapN(source[j], source[j+1]);
  435.           SwapN(aspect[j], aspect[j+1]);
  436.           SwapN(dest[j], dest[j+1]);
  437.           SwapR(&time[j], &time[j+1]);
  438.           SwapN(sign[j], sign[j+1]);
  439.           SwapN(isret[j], isret[j+1]);
  440.           j--;
  441.         }
  442.       }
  443.  
  444.       /* Now loop through list and display all the transits. */
  445.  
  446.       for (i = 0; i < occurcount; i++) {
  447.         s1 = (_int)time[i]/24/60;
  448.         s3 = (_int)time[i]-s1*24*60;
  449.         s2 = s3/60;
  450.         s3 = s3-s2*60;
  451.         SetCI(ciSave, Mon2, s1+1, Yea2, DegToDec((real)
  452.           ((_int)time[i]-s1*24*60) / 60.0), Dst2, Zon2, Lon2, Lat2);
  453.         sprintf(sz, "%s %s ",
  454.           SzDate(Mon2, s1+1, Yea2, fFalse), SzTime(s2, s3)); PrintSz(sz);
  455.         PrintAspect(source[i], sign[i], isret[i], aspect[i],
  456.           dest[i], SFromZ(planet3[dest[i]]), (int)RSgn(ret3[dest[i]]),
  457.           (char)(fProg ? 'u' : 't'));
  458.  
  459.         /* Check for a Solar, Lunar, or any other return. */
  460.  
  461.         if (aspect[i] == aCon && source[i] == dest[i]) {
  462.           AnsiColor(kWhite);
  463.           sprintf(sz, " (%s Return)", source[i] == oSun ? "Solar" :
  464.             (source[i] == oMoo ? "Lunar" : szObjName[source[i]]));
  465.           PrintSz(sz);
  466.         }
  467.         PrintL();
  468. #ifdef INTERPRET
  469.         if (us.fInterpret)
  470.           InterpretTransit(source[i], aspect[i], dest[i]);
  471. #endif
  472.         AnsiColor(kDefault);
  473.       }
  474.     }
  475.   }
  476. }
  477.  
  478.  
  479. /* Display a list of planetary rising times relative to the local horizon */
  480. /* for the day indicated in the chart information, as specified with the  */
  481. /* -Zd switch. For the day, the time each planet rises (transits horizon  */
  482. /* in East half of sky), sets (transits horizon in West), reaches its     */
  483. /* zenith point (transits meridian in South half of sky), and nadirs      */
  484. /* transits meridian in North), is displayed.                             */
  485.  
  486. void ChartInDayHorizon()
  487. {
  488.   char sz[cchSzDef];
  489.   int source[MAXINDAY], type[MAXINDAY], sign[MAXINDAY],
  490.     fRet[MAXINDAY], occurcount, division, div, s1, s2, i, j, fT;
  491.   real time[MAXINDAY], rgalt1[objMax], rgalt2[objMax], azialt[MAXINDAY],
  492.     azi1, azi2, alt1, alt2, lon, lat, mc1, mc2, xA, yA, xV, yV, d, k;
  493.  
  494.   fT = us.fSiderial; us.fSiderial = fFalse;
  495.   lon = RFromD(Mod(Lon)); lat = RFromD(Lat);
  496.   division = us.nDivision * 4;
  497.   occurcount = 0;
  498.  
  499.   ciCore = ciMain; ciCore.tim = 0.0;
  500.   CastChart(fTrue);
  501.   mc2 = RFromD(planet[oMC]); k = RFromD(planetalt[oMC]);
  502.   EclToEqu(&mc2, &k);
  503.   for (i = 1; i <= cSign; i++) {
  504.     cp2.cusp[i] = house[i];
  505.     cp2.house[i] = inhouse[i];
  506.   }
  507.   for (i = 1; i <= cObj; i++) {
  508.     cp2.obj[i] = planet[i];
  509.     rgalt2[i] = planetalt[i];
  510.     cp2.dir[i] = ret[i];
  511.   }
  512.  
  513.   /* Loop thorough the day, dividing it into a certain number of segments. */
  514.   /* For each segment we get the planet positions at its endpoints.        */
  515.  
  516.   for (div = 1; div <= division; div++) {
  517.     ciCore = ciMain; ciCore.tim = DegToDec(24.0*(real)div/(real)division);
  518.     CastChart(fTrue);
  519.     mc1 = mc2;
  520.     mc2 = RFromD(planet[oMC]); k = RFromD(planetalt[oMC]);
  521.     EclToEqu(&mc2, &k);
  522.     for (i = 1; i <= cSign; i++) {
  523.       cp1.cusp[i] = cp2.cusp[i]; cp1.house[i] = cp2.house[i];
  524.       cp2.cusp[i] = house[i];    cp2.house[i] = inhouse[i];
  525.     }
  526.     for (i = 1; i <= cObj; i++) {
  527.       cp1.obj[i] = cp2.obj[i]; cp2.obj[i] = planet[i];
  528.       rgalt1[i] = rgalt2[i]; rgalt2[i] = planetalt[i];
  529.       cp1.dir[i] = cp2.dir[i]; cp2.dir[i] = ret[i];
  530.     }
  531.  
  532.     /* For our segment, check to see if each planet during it rises, sets, */
  533.     /* reaches its zenith, or reaches its nadir.                           */
  534.  
  535.     for (i = 1; i <= cObj; i++) if (!ignore[i] && FThing(i)) {
  536.       EclToHorizon(&azi1, &alt1, cp1.obj[i], rgalt1[i], lon, lat, mc1);
  537.       EclToHorizon(&azi2, &alt2, cp2.obj[i], rgalt2[i], lon, lat, mc2);
  538.       j = 0;
  539.  
  540.       /* Check for transits to the horizon. */
  541.       if ((alt1 > 0.0) != (alt2 > 0.0)) {
  542.         d = RAbs(alt1)/(RAbs(alt1)+RAbs(alt2));
  543.         k = Mod(azi1 + d*MinDifference(azi1, azi2));
  544.         j = 1 + 2*(MinDistance(k, rDegHalf) < rDegQuad);
  545.  
  546.       /* Check for transits to the meridian. */
  547.       } else if (RSgn(MinDifference(azi1, rDegQuad)) !=
  548.         RSgn(MinDifference(azi2, rDegQuad))) {
  549.         j = 2 + 2*(MinDistance(azi1, rDegQuad) < rDegQuad);
  550.         d = RAbs(azi1 - (j > 2 ? rDegQuad : 270.0))/MinDistance(azi1, azi2);
  551.         k = alt1 + d*(alt2-alt1);
  552.       }
  553.       if (j && occurcount < MAXINDAY) {
  554.         source[occurcount] = i;
  555.         type[occurcount] = j;
  556.         time[occurcount] = 24.0*((real)(div-1)+d)/(real)division*60.0;
  557.         sign[occurcount] = (int)Mod(cp1.obj[i] +
  558.           d*MinDifference(cp1.obj[i], cp2.obj[i]))/30 + 1;
  559.         fRet[occurcount] = (int)RSgn(cp1.dir[i]) + (int)RSgn(cp2.dir[i]);
  560.         azialt[occurcount] = k;
  561.         ciSave = ciMain;
  562.         ciSave.tim = DegToDec(time[occurcount] / 60.0);
  563.         occurcount++;
  564.       }
  565.     }
  566.   }
  567.  
  568.   /* Sort each event in order of time when it happens during the day. */
  569.  
  570.   for (i = 1; i < occurcount; i++) {
  571.     j = i-1;
  572.     while (j >= 0 && time[j] > time[j+1]) {
  573.       SwapN(source[j], source[j+1]);
  574.       SwapN(type[j], type[j+1]);
  575.       SwapR(&time[j], &time[j+1]);
  576.       SwapN(sign[j], sign[j+1]);
  577.       SwapN(fRet[j], fRet[j+1]);
  578.       SwapR(&azialt[j], &azialt[j+1]);
  579.       j--;
  580.     }
  581.   }
  582.  
  583.   /* Finally display the list showing each event and when it occurs. */
  584.  
  585.   for (i = 0; i < occurcount; i++) {
  586.     ciSave = ciMain;
  587.     ciSave.tim = DegToDec(time[i] / 60.0);
  588.     j = DayOfWeek(Mon, Day, Yea);
  589.     AnsiColor(kRainbowA[j + 1]);
  590.     sprintf(sz, "(%c%c%c) ", chDay3(j)); PrintSz(sz);
  591.     AnsiColor(kDefault);
  592.     s1 = (int)time[i]/60;
  593.     s2 = (int)time[i]-s1*60;
  594.     sprintf(sz, "%s %s ", SzDate(Mon, Day, Yea, fFalse), SzTime(s1, s2));
  595.     PrintSz(sz);
  596.     AnsiColor(kObjA[source[i]]);
  597.     sprintf(sz, "%7.7s ", szObjName[source[i]]); PrintSz(sz);
  598.     AnsiColor(kSignA(sign[i]));
  599.     sprintf(sz, "%c%c%c%c%c ",
  600.       fRet[i] > 0 ? '(' : (fRet[i] < 0 ? '[' : '<'), chSig3(sign[i]),
  601.       fRet[i] > 0 ? ')' : (fRet[i] < 0 ? ']' : '>')); PrintSz(sz);
  602.     AnsiColor(kElemA[type[i]-1]);
  603.     if (type[i] == 1)
  604.       PrintSz("rises  ");
  605.     else if (type[i] == 2)
  606.       PrintSz("zeniths");
  607.     else if (type[i] == 3)
  608.       PrintSz("sets   ");
  609.     else
  610.       PrintSz("nadirs ");
  611.     AnsiColor(kDefault);
  612.     PrintSz(" at ");
  613.     if (type[i] & 1) {
  614.       j = (int)(RFract(azialt[i])*60.0);
  615.       sprintf(sz, "%3d%c%02d'", (int)azialt[i], chDeg1, j); PrintSz(sz);
  616.  
  617.       /* For rising and setting events, we'll also display a direction  */
  618.       /* vector to make the 360 degree azimuth value thought of easier. */
  619.  
  620.       xA = RCosD(azialt[i]); yA = RSinD(azialt[i]);
  621.       if (RAbs(xA) < RAbs(yA)) {
  622.         xV = RAbs(xA / yA); yV = 1.0;
  623.       } else {
  624.         yV = RAbs(yA / xA); xV = 1.0;
  625.       }
  626.       sprintf(sz, " (%.2f%c %.2f%c)",
  627.         yV, yA < 0.0 ? 's' : 'n', xV, xA > 0.0 ? 'e' : 'w'); PrintSz(sz);
  628.     } else
  629.       PrintAltitude(azialt[i]);
  630.     PrintL();
  631.   }
  632.  
  633.   /* Recompute original chart placements as we've overwritten them. */
  634.  
  635.   ciCore = ciMain;
  636.   us.fSiderial = fT;
  637.   CastChart(fTrue);
  638. }
  639.  
  640.  
  641. /* Print out an ephemeris - the positions of the planets (at the time in the */
  642. /* current chart) each day during a specified month, as done with the -E     */
  643. /* switch. Display the ephemeris for the whole year if -Ey is in effect.     */
  644.  
  645. void ChartEphemeris()
  646. {
  647.   char sz[cchSzDef];
  648.   int yea, yea1, yea2, mon, mon1, mon2, daysiz, i, j, s, d, m;
  649.  
  650.   /* If -Ey is in effect, then loop through all months in the whole year. */
  651.  
  652.   if (us.nEphemYears) {
  653.     yea1 = Yea; yea2 = yea1 + us.nEphemYears - 1; mon1 = 1; mon2 = 12;
  654.   } else {
  655.     yea1 = yea2 = Yea; mon1 = mon2 = Mon;
  656.   }
  657.  
  658.   /* Loop through the year or years in question. */
  659.  
  660.   for (yea = yea1; yea <= yea2; yea++)
  661.  
  662.   /* Loop through the month or months in question, printing each ephemeris. */
  663.  
  664.   for (mon = mon1; mon <= mon2; mon++) {
  665.     daysiz = DayInMonth(mon, yea);
  666.     PrintSz(us.fEuroDate ? "Dy/Mo/Yr" : "Mo/Dy/Yr");
  667.     for (j = 1; j <= cObj; j++) {
  668.       if (!ignore[j] && FThing(j)) {
  669.         sprintf(sz, "  %s%c%c%c%c%s ", is.fSeconds ? "  " : "", chObj3(j),
  670.           szObjName[j][3] != 0 ? szObjName[j][3] : ' ',
  671.           is.fSeconds ? "   " : ""); PrintSz(sz);
  672.       }
  673.     }
  674.     PrintL();
  675.     for (i = 1; i <= daysiz; i = AddDay(mon, i, yea, 1)) {
  676.  
  677.       /* Loop through each day in the month, casting a chart for that day. */
  678.  
  679.       SetCI(ciCore, mon, i, yea, Tim, Dst, Zon, Lon, Lat);
  680.       CastChart(fTrue);
  681.       PrintSz(SzDate(mon, i, yea, -1));
  682.       PrintCh(' ');
  683.       for (j = 1; j <= cObj; j++)
  684.         if (!ignore[j] && FThing(j)) {
  685.           if (is.fSeconds)
  686.             PrintZodiac(planet[j]);
  687.           else {
  688.             AnsiColor(kObjA[j]);
  689.             s = SFromZ(planet[j]);
  690.             d = (int)planet[j] - (s-1)*30;
  691.             m = (int)(RFract(planet[j])*60.0);
  692.             sprintf(sz, "%2d%s%02d", d, szSignAbbrev[s], m); PrintSz(sz);
  693.           }
  694.           PrintCh((char)(ret[j] >= 0.0 ? ' ' : '.'));
  695.         }
  696.       PrintL();
  697.       AnsiColor(kDefault);
  698.     }
  699.     if (mon < mon2 || yea < yea2)
  700.       PrintL();
  701.   }
  702. }
  703.  
  704. /* charts3.c */
  705.